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In  this  paper  we  report  on  the  development  and  demonstration  of  physically  consistent  three- 
dimensional  models  for  Lithium  Ion  Battery  (LIB)  cells.  The  discharge  behavior  of  a  LIB  is  a  multi¬ 
physics  and  multiscale  problem  that  is  simulated  using  coupled  models  for  thermal,  electrical,  and 
electrochemical  phenomena.  The  individual  physics  models  and  software  are  integrated  into  a  new  open 
computational  framework  for  battery  simulations  which  was  designed  to  support  a  variety  of  modeling 
formulations  and  computer  codes.  Several  cell  configurations  (unrolled  cell,  unrolled  cell  with  current 
collectors,  large  capacity  pouch  cell,  and  cylindrical  cell)  that  show  the  importance  of  coupled  simula¬ 
tions  are  simulated  using  this  approach  and  discussed.  A  validation  study  is  presented  for  the  pouch  cell 
discharged  under  high  rates  to  demonstrate  the  accuracy  of  the  proposed  modeling  framework. 

©  2013  Elsevier  B.V.  All  rights  reserved. 


1.  Introduction 

High  energy  density  lithium-ion  batteries  are  the  primary  po¬ 
wer  source  choice  for  electric  vehicles  and  large-scale  military  and 
space  applications  [1].  During  various  abuse  conditions,  over¬ 
heating  and  uneven  cooling  of  the  battery  and  battery  pack  can 
result  in  rapid  degradation  and  in  extreme  cases,  uncontrolled 
electro-thermo-chemical  reaction.  High-performance  battery 
packs  or  cells  are  designed  to  best  perform  within  narrow  tem¬ 
perature  ranges.  Mathematical  modeling  of  battery  thermal 
behavior  and  cooling  strategies  has  proven  to  be  an  efficient  and 
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cost  effective  design  tool  for  predicting  life  cycle  issues  and  failure 
modes  2,3].  In  order  to  make  accurate  life  predictions,  to  improve 
performance  and  extend  durability  of  the  battery,  the  thermal 
variations  within  the  electrode  material  must  be  taken  into 
consideration.  Approaches  based  on  homogenization  of  the  cell 
with  averaged  properties  may  not  be  appropriate  when  the  thermal 
runaway  and  (or)  electrical  shorts  need  to  be  modeled  due  to  lack  of 
proper  resolution  of  current  collectors. 

Several  authors  [4-6  have  recently  proposed  coupling  between 
the  thermal  models  and  representative  electric  circuit  models  for 
battery  module.  The  coupling  is  implemented  in  two  steps.  First,  an 
equivalent  electric  circuit  is  constructed  to  model  the  electric  cur¬ 
rents  and  overpotentials  in  the  module.  Then,  the  electrical  solution 
is  used  to  calculate  the  sources  for  the  thermal  equation.  In  the 
thermal  model,  the  entire  cell  is  homogenized  based  on  volume 
averaged  physical  properties  (e.g.,  thermal  conductivity,  thermal 
capacity,  density)of  the  electrodes  and  current  collectors.  While 
being  three-dimensional,  such  approach  loses  the  spatial  resolution 
of  the  currents  in  the  collectors  due  to  such  circuit  homogenization. 
In  Ref.  [7],  authors  propose  a  new  3D  multiphysics  model  for  Li-ion 
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battery,  where  the  electric  currents  in  collectors  are  modeled  using 
two  variables  (one  for  positive  current  collector  and  the  other  one 
for  negative).  The  currents  are  coupled  with  the  electrochemical 
model  via  surface  current  density  and  the  anisotropies  in  the 
connecting  tabs  are  modeled  using  geometric  relations. 

In  order  to  understand  and  predict  the  complex  interplay  be¬ 
tween  various  coupled  physics  phenomena  in  Li-ion  batteries,  many 
electrochemical,  thermal,  and  life  models  were  developed  [8-14]. 
Several  codes  [15-18]  have  been  developed  to  bring  these  models 
under  a  single  framework.  However  all  these  codes  require  the  user 
to  either  rewrite  or  develop  the  model  within  the  framework  or  to 
use  the  models  that  pre-exist  in  the  code.  Currently,  there  are  no 
open-source  computational  platforms  that  can  integrate  these  in¬ 
dividual  models  beyond  their  original  scope  and  application.  Given 
the  complex  requirements  for  batteries  for  transportation  applica¬ 
tions,  a  predictive  simulation  capability  which  can  guide  design  by 
considering  performance  and  safety  implications  of  different 
chemistry  and  materials  choices  is  needed.  We  have  developed  an 
open  computational  framework  that  integrates  the  existing  models 
that  are  independently  developed  and  can  also  accommodate  new 
models  for  simulating  battery  performance.  The  framework  and  the 
underlying  models  span  across  the  cells,  modules,  battery  pack,  etc. 

We  are  proposing  an  approach  where  electrochemistry,  thermal 
and  electrical  transport  are  solved  within  a  single  simulation 
framework  using  three-dimensional  representation.  Multi-material 
interfaces  are  kept  intact  without  homogenization  assumption.  A 
schematic  for  the  cell  sandwich  model  is  shown  in  Fig.  1.  Domains  S 
and  P  represent  cell  sandwich  (anode,  separator,  and  cathode)  and 
pouch  material;  PC  and  NC  represent  positive  and  negative  current 
collectors.  In  Section  2,  we  describe  the  mathematical  models  for 
various  physics  and  software  components  that  are  used  in  the 
framework  for  coupled  simulation.  In  Section  3,  we  present  the 
computing  infrastructure  and  the  underlying  numerical  methods 
for  the  framework.  In  Section  4,  several  battery  simulation  exam¬ 
ples  are  presented  to  demonstrate  the  capability  and  accuracy  of 
the  approach. 

2.  Mathematical  models 

One  of  our  main  goals  is  to  develop  a  modular  framework  so  that 
individual  components  that  correspond  to  particular  physical 
phenomena  can  be  flexibly  and  consistently  interchanged  with 
similar  components.  In  this  section,  we  describe  the  models  for 
underlying  physics  phenomena,  their  mathematical  formulations 
and  computational  implementations.  The  coupling  between  the 
models  is  discussed  in  Section  3. 

2.1.  Electrochemical  models 

Two  models  for  electrochemical  behavior  were  chosen  in  the 
present  work.  The  first  model  takes  experimentally  measured  po¬ 
larization  characteristics  of  the  electrochemical  cell  as  an  input. 
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Fig.  2.  Determination  of  the  OCP  and  polarization  parameter  from  the  discharge  data. 

Linear  dependency  of  discharge  cell  voltage  on  current  density  is 
given  by  Eq.  (1)  from  Newman,  Tiedemann  19]  and  Gu  [20]  for  an 
electrochemical  cell.  This  model  is  as  also  referred  as  NTG  (New¬ 
man,  Tiedemann,  and  Gu)  model.  The  equations  for  the  model  are: 

(J  =  Y  (Vp-Vn-U) 

u-|y 

M 

y=E  bid1 

i  =  0 

Variables  Y  and  U  represent  effective  conductance  and  open  cir¬ 
cuit  potential  of  the  cell  as  a  function  of  depth  of  discharge  ( 6 ),  and  J 
is  the  current  density  transferred  from  the  negative  electrode  to  the 
positive  electrode.  Constants  a*  and  bi  are  fitting  parameters  that  are 
determined  from  cell  discharge  curves  for  a  number  of  different  C- 
rates.  When  the  cell  potentials  are  plotted  against  the  current  den¬ 
sity  coordinates,  the  intercepts  of  these  curves  at  zero  current 
density  represent  the  values  of  open  circuit  potential  U  and  the 
slopes  represent  the  reciprocal  of  the  cell  polarization  parameter  Y. 
This  is  shown  in  Fig.  2,  where  the  experimental  data  from  the  pouch 
cell  discussed  further  in  Section  4.3  is  presented  to  aid  the 
description  of  the  method.  Further  details  of  the  parameter  cali¬ 
bration  procedure  can  be  found  in  Refs.  [19,20].  The  model  Eq.  (1) 
describes  the  overall  characteristics  of  the  cell  without  resolving 
local  kinetics  of  electrochemical  reactions  or  the  internal  concen¬ 
tration  and  potential  gradients  within  the  cell-sandwich.  It  is  a 
linearized  approximation  of  the  full  polarization  curve.  The  NTG 
model  has  been  widely  used  in  cell  and  pack-level  thermal  modeling 
[3,21,22]  due  to  simplicity  of  its  implementation,  provided  the  cell 
discharge  curves  are  known  for  sufficient  number  of  discharge  rates. 

The  second  implemented  model  for  electrochemistry  is  based 
on  the  porous  electrode  theory  [13,23-25].  The  corresponding 


Fig.  1.  Schematic  of  a  Li-ion  cell. 
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component  within  the  framework  is  called  DualFoil.  The  model 
employs  one-dimensional  description  of  the  Li-ion  transport  from 
the  negative  electrode  and  into  the  positive  electrode.  Transport 
through  the  electrolyte  is  modeled  by  using  the  concentrated  so¬ 
lution  theory  resulting  in  the  equation  for  the  lithium  concentra¬ 
tion  (ce)  along  the  transversal  direction  z  within  the  cell  sandwich 


where  e  is  the  volume  fraction  of  electrolyte  (assumed  to  be  time- 
independent)  and  is  equal  to  1  in  the  separator  region.  Dfff  is  the 
effective  diffusivity  of  lithium  ion  in  the  electrolyte,  is  the 
cation  transference  number,  and  jn  is  the  flux  of  Li+  across  the 
interface  between  the  electrolyte  and  active  electrode  material. 
Obviously  jn  =  0  in  the  separator  region.  The  lithium  flux  is  con¬ 
nected  with  the  cell  overpotential  via  the  Butler— Volmer  kinetics 
equation 


2.2.  Thermal  transport  model 

The  transient  three-dimensional  heat  conduction  equation 
derived  from  the  thermal  energy  conservation  is: 

pCp^-V(kvr)  =  q  (7) 

Where  p  is  the  density,  Cp  is  the  specific  heat  capacity,  k  =  {  kx,  ky, 
kz)  is  the  thermal  conductivity,  and  T  is  local  temperature.  A 
general  form  of  the  heat  generation  term  q  in  a  battery  has  been 
derived  by  Bernardi  et  al.  based  on  the  energy  balance  within  the 
cell  [27].  The  change  in  temperature  within  a  battery  cell  is 
assumed  to  be  caused  by  electrochemical  reactions,  changes  in  the 
heat  capacity  of  the  system,  phase  changes,  mixing,  electrical 
work,  and  heat  transfer  with  the  surroundings.  Typically  the  heat 
of  mixing  and  phase  changes  is  ignored  for  the  electrochemical 
systems  [3,25,27-29],  and  the  simplified  form  of  the  heat  gener¬ 
ation  with  addition  of  ohmic  heating  in  the  electrodes  can  be 
written  as 


-  exp 


(3) 


where  i0  is  the  exchange  current  density,  aa,  ac  are  the  anodic  and 
cathodic  transfer  coefficients  respectively,  and  p  is  the  over¬ 
potential  (pj  =  0s  -  Qe  -  ^Q’.ocp))-  represent  the  local  potential 
of  the  electrolyte  and  the  electrode  material  respectively,  and  Ujtocp 
is  the  open  circuit  potential  for  anode  and  cathode  as  a  function  of 
lithium  concentration  in  solid  phase  (cs). 

The  concentration  of  lithium  in  the  solid  phase  (cs)  is  deter¬ 
mined  by  the  diffusion  equation 


9cs 
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]_d_ 

r2  dr 


Dfr 2^5 

s  gr 


=  0 


(4) 


with  r  being  the  radial  coordinate  of  a  representative  spherical 
active  material  particle  and  representing  the  effective  diffu¬ 
sivity.  The  above  transport  equations  are  complemented  by  the 
equations  corresponding  to  the  conservation  of  charge  in  electro¬ 
lyte  Eq.  (5)  and  in  active  material  Eq.  (6) 


<7  =  Y)  -<t>e-  Uj)  +  aihT  -Qf  +  ffeffV<fe  •  V0S 
j  j 

+  keffV^  •  V^e  +  (l  -  t°  )  (l  +  Vln  Ce  •  V</>e 

(8) 


where  the  first  term  is  the  irreversible  heat  energy  generated  due 
to  cell  polarization,  and  the  second  term  represents  the  revers¬ 
ible  entropy  change  due  to  a  reaction  j.  Quite  often  this  latter 
entropic  term  is  set  to  zero,  either  due  to  the  lack  of  data 
describing  the  open  circuit  potential  as  a  function  of  temperature 
[28  or  based  on  the  assumption  of  negligible  influence  on  the 
total  heat  generation  [30].  The  rest  of  the  terms  in  Eq.  (8) 
represent  the  effect  of  ohmic  heating  in  solid  and  electrolyte 
phases.  The  latter  is  typically  ignored  [21,22,28]  so  that  the  Eq. 
(8)  simplifies  to 


Q  —  ajij  (0s  0e 
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Eq.  (9)  is  complemented  by  the  sources  due  to  ohmic  heating  in 

(5)  current  collectors  which  are  obtained  from  the  electrical  compo¬ 
nent  (discussed  in  the  next  section). 

When  the  NTG  model  is  used  within  the  electrochemical 
component,  the  heat  source  is  represented  as  an  integral  quantity 

(6)  of  the  source  in  Eq.  (8)  in  form 


where  is  the  molar  activity  coefficient  of  the  salt  in  the  elec¬ 
trolyte.  keff  and  a  represent  respectively  the  ionic  conductivity  of 
the  electrolyte  and  electronic  conductivity  of  the  solid  active  ma¬ 
terial.  The  second  term  in  Eq.  (5)  represents  contribution  of  po¬ 
larization  and  obviously  Eq.  (5)  reduces  to  Kirchhoffs  law  Eq.  (6) 
for  electrolytes  with  unity  transference  number.  Influence  of  the 
transference  number  of  electrolyte  is  discussed  in  detail  in  Ref. 
[26]. 

As  can  be  seen,  the  entire  electrochemical  system  can  be  cast 
into  four  conservation  equations.  Transport  in  electrolyte  and 
electrode  is  given  by  Eqs.  (2)— (4),  and  the  conservation  of  charge  is 
given  by  Eqs.  (5)  and  (6).  In  this  model  it  is  assumed  that  the 
porosity  of  the  electrode  remains  constant  and  that  the  volumetric 
strains  resulting  from  lithium  insertion  are  negligible.  The  details  of 
the  model  derivation  can  be  found  in  Refs.  [13,23-25]. 


Q  = 


1 

h 


(10) 


where  rj  is  the  cell  overpotential  ( Vp  -  Vn  -  U),  J  is  the  current 
density  passing  through  the  cell  (applied  current  normalized  with 
respect  to  the  total  cell  surface  area),  and  h  is  the  cell  thickness. 
Such  integral  representation  is  because  NTG  model  is  based  on 
homogenization  across  the  entire  cell  and  thus  cannot  resolve  the 
internal  profiles  of  the  lithium  concentration  or  electrolyte/elec¬ 
trode  potentials. 

The  components  of  effective  thermal  conductivity  in  Eq.  (7)  are 
calculated  taking  into  account  the  anisotropy  of  the  thermal 
properties  of  the  cell  sandwich.  With  z-axis  along  the  thickness  of 
the  cell  sandwich,  the  transversal  and  longitudinal  components  of 
the  effective  thermal  conductivity  are  calculated  as 
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'kz  = 


kx  = 


(11) 


where  ft,  and  /q  are  the  individual  component  thickness  and  ther¬ 
mal  conductivities. 


2.3.  Electrical  transport  model 

Electrical  component  within  the  framework  is  based  on  the 
Laplace’s  equation  (equation  for  conservation  of  charge:  Kirchh- 
off  s  law) 

V(<7  W)  =  0  (12) 

The  effective  cell  sandwich  conductivity  (provided  by  the  elec¬ 
trochemistry  model  in  domain  S  of  Fig.  1 ),  and  current  collector 
conductivity  (in  domains  NC  and  PC)  is  used  when  solving  Laplace’s 
equation  on  a  single  3D  domain  with  multi-material  interfaces  to 
maintain  a  consistent  potential  drop  in  the  electrical  simulation. 
For  discharge  behavior,  a  constant  current  is  imposed  as  boundary 
condition  on  the  edge  of  the  current  collector  tabs. 

Heat  transfer  equations  and  Laplace’s  equation  for  electronic 
transport  are  solved  using  the  Amperes,  a  multiphysics  3D  Finite 
Element  (FE)  code  developed  at  ORNL.  Given  the  different 
spatiotemporal  scale  of  electrochemical  and  thermal  transport, 
we  undertake  the  simulations  at  respective  scales  using  operator 
splitting  method  so  that  we  have  an  efficient  computational 
scheme  to  solve  the  coupled  problem.  The  electrochemistry 
model  requires  a  more  refined  discretization  along  the  thickness 
of  the  cell  compared  to  the  thermal  transport  where  a  coarser 
discretization  is  sufficient. 

The  information  between  the  different  discretizations  and 
models  is  exchanged  by  decomposing  the  cell  domain  into  multiple 
zones  along  the  various  coordinate  axes  to  account  for  the  effect  of 
temperature  gradients  on  the  electrochemical  model.  In  each  zone, 
an  instance  of  electrochemical  simulation  is  launched.  The  elec¬ 
trochemical  solution  provides  the  heat  source  for  each  of  the  zones. 
The  Amperes  component  updates  the  temperatures  for  the  each 
of  the  zones  for  the  subsequent  run  of  the  electrochemical 
component. 


3.  Computational  framework 

The  Open  Architecture  Software  (OAS)  framework  for  battery 
simulations  was  developed  under  the  US  Department  of  Energy 
research  program  Computer  Aided  Engineering  for  Batteries 
(CAEBAT).  The  framework  is  designed  for  transient  simulations 
with  relatively  loose  coupling  between  the  various  physics  com¬ 
ponents.  Each  simulation  component  in  the  OAS  framework  is 
equipped  with  coupling  interfaces  that  are  implemented  as  python 
wrappers  around  standalone  code  executables.  These  interfaces 
convert  native  data  from  these  executables  to  the  data  format  used 
by  the  OAS  framework.  This  enables  the  executables  to  remain 
unchanged  and  yet  fully  interact  with  the  framework  and  other 
components.  The  OAS  framework,  the  physics  components,  in¬ 
terfaces  and  the  support  tools,  constitute  the  Virtual  Integrated 
Battery  Environment  (VIBE). 

With  this  framework,  scientists  and  designers  can  utilize  de¬ 
cades  of  effort  that  various  communities  have  put  into  the  devel¬ 
opment,  verification,  and  validation  of  modeling  tools  for  the 
individual  physics  domains,  and  focus  on  the  new  issues  associated 
with  coupling  and  multiphysics  interaction.  In  many  cases,  the  data 
exchanged  by  components  is  modest,  and  is  easily  managed 
through  a  battery  state.  The  battery  state  is  the  minimal  digital 
description  of  the  battery  in  space  and  time  such  that  each  simu¬ 
lation  component  can  apply  their  respective  physics  models  and 
advance  in  time  from  one  state  point  to  the  next. 

The  architecture  of  the  OAS  framework  is  shown  in  Fig.  3.  The  OAS 
framework  is  derived  from  the  Integrated  Plasma  Simulator  (IPS) 
developed  by  the  Center  for  Simulation  of  RF  Wave  Interactions  with 
Magneto-hydrodynamics  (SWIM)  [31].  The  OAS  is  designed  pri¬ 
marily  for  use  in  a  batch  processing  environment,  with  a  batch  job 
typically  comprising  a  single  invocation  of  the  framework,  calling  the 
individual  physics  codes  underlying  the  components  many  times  as 
the  simulation  progresses.  OAS  simulations  are  typically  orches¬ 
trated  by  a  “driver”  component,  though  simulation-controlling  logic 
can  also  be  built  into  individual  components.  The  OAS  framework 
provides  a  set  of  services  that  are  used  by  components  within  a 
simulation  to  facilitate  the  coupling  process.  The  “framework  ser¬ 
vices”  communicates  with  the  various  physics  based  components 
through  “component  adapters.”  These  components  update  the 
“battery  state”  through  “state  adapters.”  The  detailed  description  of 
the  functionality  and  services  provided  by  the  OAS  framework  can  be 
found  in  Refs.  [32,33]. 


Fig.  3.  Schematic  of  OAS  framework. 
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The  OAS  framework  supports  multiple  levels  of  concurrency, 
enabling  efficient  and  flexible  utilization  of  heterogeneous 
computational  resources  [32].  Individual  physics  based  compo¬ 
nents  can  be  serial  or  parallel,  multiple  instances  of  the  component 
tasks  can  be  executed  concurrently  in  a  simulation,  and  multiple 
simulations  can  also  be  run  simultaneously.  To  demonstrate  these 
features,  we  have  implemented  components  into  VIBE  that  support 
various  levels  of  parallelism  for  modeling  electrochemistry  and 
mass,  electron,  and  heat  transport.  The  schematic  of  these  com¬ 
ponents  interacting  with  the  battery  state  is  illustrated  in  Fig.  4. 

Based  on  the  interaction  of  the  dependent  variables  and  the 
coupling,  various  computational  strategies  between  these  systems 
can  be  employed.  When  the  influence  between  two  states  is  pre¬ 
dominantly  in  one  direction,  a  one-way  coupled  numerical  tech¬ 
nique  suffices  to  get  plausible  solutions.  In  a  one-way  or  forward 
coupling,  a  set  of  variables  in  one  of  the  physics  components  de¬ 
pends  on  the  solution  of  the  other  components,  but  not  vice  versa. 


On  the  other  hand,  as  the  inter-dependence  between  the  variables 
becomes  strong  a  two-way  coupled  numerical  technique  must  be 
used  to  obtain  an  accurate  solution.  In  this  case,  the  state  variables 
within  the  components  strongly  depend  on  the  solutions  of  the 
other  components  and  exchange  needs  to  occur  in  both  directions. 
Picard  iterative  method  is  applicable  to  fully  coupled  physics,  for 
example  in  cases  when  the  interfaces  or  domains  are  shared  be¬ 
tween  the  two  sets  of  physics.  In  the  case  of  explicit  coupling,  the 
exchange  of  data  between  different  physics  occurs  at  the  end  of  the 
time  step,  while  the  implicit  coupling  imposes  simultaneous  solu¬ 
tion  of  all  physics.  Various  coupling  scenarios  available  for  battery 
modeling  are  described  in  Fig.  5. 

4.  Numerical  studies 

To  illustrate  the  flexibility  of  using  various  components  under 
OAS  framework  and  effect  of  choice  of  representative  physics 
models,  a  series  of  three  loosely  coupled  simulations  of  electro¬ 
chemistry,  electrical  and  thermal  transport  are  discussed  for 
several  battery  configurations.  The  developed  framework  was  used 
to  simulate  multiphysics  processes  in  typical  battery  cell  configu¬ 
rations.  The  coupling  between  different  physics  solvers  was 
implemented  as  a  two-way  loose  coupling  according  to  Fig.  5  and 
the  number  of  interaction  points  between  the  two  components  has 
been  varied  to  reach  convergence.  A  validation  study  demonstrates 
a  good  fit  between  the  experimentally  measured  temperature 
profile  and  the  simulation.  All  the  simulations  were  executed  on  8 
core  AMD  processor  workstation  with  8  GB  memory. 

4.1.  Unrolled  cell 

In  the  first  case  study,  we  simulate  the  unrolled  cell  as  described 
in  Ref.  [29].  Cell  electrochemistry  was  modeled  by  DualFoil 
component.  Two  simulation  scenarios  were  considered.  In  the  first, 
DualFoil  was  coupled  with  thermal  transport  only.  The  second 
scenario  couples  DualFoil  to  thermal  and  electrical  models. 

The  composition  of  the  cell  is  identical  to  that  described  in  Refs. 
[23,29].  Lithium  manganese  oxide  was  used  as  an  active  cathode 
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Fig.  5.  Coupling  scenarios  in  battery  modeling. 
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material  and  petroleum  derived  carbon  coke  was  used  as  anode. 
The  electrolyte  consisted  of  LiPF6  salt  in  non-aqueous  mixture  of 
ethylene  carbonate  (EC)  and  dimethyl  carbonate  (DMC)  in  a 
poly(vinyl  difluoride)-hexafluoropropylene  (PVDF-HFP)  polymer 
matrix.  The  coefficients  of  the  equation  describing  the 
concentration-dependent  conductivity  of  the  electrolyte,  as  well  as 
the  other  model  parameters  (transference  numbers,  electrode  film 
thicknesses,  etc.)  were  adopted  from  Ref.  [23].  Effective  ionic  con¬ 
ductivity  (keff)  in  Eq.  (5)  follows  the  effective  medium  approach  as 

keff  =  epke  (13) 

with  p  =  3.3  and  the  conductivity  of  electrolyte  represented  as  a 
function  of  the  lithium  concentration  [23] 

ke  =  4.1253  10-4  +  5.007  10_3ce  —  4.7212- 10-3qr 

+  1.5094  10_3c|  -  1.6018- 10_4c4  (14) 

The  effective  diffusivity  (D£tr)  in  Eq.  (2)  follows  the  same  effec- 
tive  medium  law  Eq.  (13)  with  the  same  exponent  equal  to  3.3.  The 
length  and  the  width  of  the  unrolled  cell  were  50  cm  and  2.4  cm, 
respectively  [29  ,  which  equates  to  1C  discharge  current  density  of 
17.5  Ah  m-2.  The  values  of  the  parameters  used  in  the  electro¬ 
chemical  model  are  listed  in  Table  1.  Table  2  lists  the  parameters 
used  for  thermal  modeling.  The  thermal  conductivity  values  for 
positive  and  negative  electrode  materials  listed  in  Table  2  were 
taken  from  Ref.  [34]. 

The  temperature  distribution  in  the  unrolled  cell  corresponding 
to  two  different  C-rates  is  shown  in  Fig.  3.  The  temperature  shown 
is  at  the  end  of  discharge.  The  time  step  size  (for  interaction  of  data 
between  the  components  while  internally  DualFoil  might  take 
much  smaller  time  steps)  sensitivity  was  tested  with  different 
values.  It  was  determined  that  the  difference  in  maximum  tem¬ 
perature  prediction  between  10  and  20  time  steps  was  0.3%  which 
indicates  good  convergence.  Cell  cooling  was  simulated  by  pre¬ 
scribing  temperature  of  298  K  at  the  top  surface.  The  remaining 
surfaces  are  modeled  as  perfect  insulators,  zero  heat  flux  or  adia¬ 
batic  [29],  which  is  the  default  boundary  condition  for  FEM  thermal 
formulation.  These  conditions  would  lead  to  thermal  gradient 
which  requires  that  the  cell  be  divided  into  15  zones  (segments) 
along  the  longitudinal  y-axis  (as  the  temperature  varies  in  that 
direction)  with  DualFoil  component  providing  electrochemical 
state  for  each  zone.  For  the  composition  of  the  cell  under  consid¬ 
eration  and  the  thermal  properties  of  the  components  15  zones 
distributed  along  longitudinal  y-axis  of  the  cell  was  determined  to 
be  sufficient  to  obtain  convergence  in  results  of  thermal  simulation. 


Table  1 

Electrochemical  ( DualFoil )  modeling  parameters. 


Parameter 

Value 

Thickness  of  negative  electrode  [m] 

128.0  x  10"6 

Thickness  of  positive  electrode  [m] 

190.0  x  10"6 

Thickness  of  negative  current  collector  [m] 

10.0  x  10"6 

Thickness  of  positive  current  collector  [m] 

15.0  x  10“6 

Initial  salt  concentration  [mol  ur3] 

2000 

Diffusion  coefficient  in  LiMn204  [m2  s-1] 

1.0  x  10-13 

Diffusion  coefficient  in  LiC6  [m2  s-1] 

3.9  x  lO"14 

Volume  fraction  of  electrolyte  in  positive  electrode 

0.444 

Volume  fraction  of  electrolyte  in  negative  electrode 

0.357 

Radius  of  particles  in  positive  electrode  [m] 

8.5  x  10"6 

Radius  of  particles  in  negative  electrode  [m] 

12.5  x  10“6 

Conductivity  of  positive  electrode  [S  m-1] 

3.8 

Conductivity  of  negative  electrode  [S  m-1] 

100.0 

Cathodic  reaction  rate  constant 

3.0  x  10"11 

Anodic  reaction  rate  constant 

1.0  x  10"5 

Cation  transference  number 

0.2 

Table  2 

Parameters  used  for  thermal  modeling. 


Component 

Density 

(kg  m-3) 

Heat  capacity 
(J  kg-1 1C1) 

Thermal  conductivity 
(W  nr1  K 1) 

Positive  current 

2700 

900 

238 

collector 

Positive  electrode 

1500 

700 

5 

Negative  current 

8960 

385 

398 

collector 

Negative  electrode 

2500 

700 

5 

Separator 

1200 

700 

1 

As  can  be  seen,  the  discharge  at  1C  (Fig.  6(a))  and  2C  (Fig.  6(b)) 
causes  significant  temperature  gradients,  which  are  greater  at  the 
higher  discharge  rate.  As  much  as  70  °C  maximum  temperature  can 
be  achieved  at  the  end  of  discharge  at  2C. 

The  lithium  concentration  profiles  in  electrolyte  for  various 
discharge  rates  across  the  cell  sandwich  are  shown  in  Fig.  7(a). 
These  are  obtained  at  the  end  of  discharge  when  the  cell  potential 
has  dropped  to  2  V  and  one  can  see  severe  gradients  in  the  Li+ 
concentration  due  to  slow  diffusion  across  the  thick  cathode.  Sig¬ 
nificant  gradients  in  electrode  composition  Fig.  7(b)  develop  with 
increasing  applied  discharge  rates  due  to  rate  of  lithium  insertion 
into  manganese  oxide  cathode  and  compounded  by  the  gradients 
in  Li+  concentration  in  the  electrolyte.  The  difference  between  the 
maximum  and  minimum  Li  content  in  LiMn204  through  the  elec¬ 
trode  thickness  is  as  much  as  4.6%  in  case  of  4C  discharge  rate. 

Next,  the  current  collectors  and  the  electrical  model  component 
Eq.  (12)  were  added  to  the  above  system  to  model  the  NC  and  PC 
domains.  The  electrical  component  provides  the  electrical  current 
density  solution  within  the  cell  sandwich  thus  providing  the  heat 
sources  for  the  thermal  component  within  the  NC  and  PC  domains. 
Fig.  8  depicts  the  temperature  distribution  of  this  unrolled  cell  with 
current  collectors.  It  should  be  mentioned,  that  while  inclusion  of 
electrical  model  and  current  collectors  results  in  additional  heat 
sources  arising  from  the  ohmic  heating,  the  magnitude  of  such 
heating  is  much  smaller  than  the  sources  coming  from  polarization 
Eq.  (9).  As  can  be  seen  from  comparison  of  Figs.  6  and  8,  addition  of 
current  collectors  results  in  overall  lower  temperatures  within  the 
cell,  due  to  the  heat  dissipation  through  metal  current  collectors. 
Compared  to  the  results  in  Fig.  6  the  maximum  temperature 
decreased  by  3%  and  the  temperature  gradient  within  the  cell 
longitudinal  direction  decreased  as  well.  The  total  time  taken  to 
complete  a  single  simulation  is  around  3  min  with  a  finite  element 


Fig.  6.  Temperature  distribution  (in  K)  in  unrolled  cell  at  the  end  of  discharge  at:  (a) 
1C,  (b)  2C  (Not  to  scale). 
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Fig.  7.  Effect  of  discharge  rate  on  (a)  lithium  concentration  in  electrolyte,  (b)  electrode 
composition. 

mesh  of  3000  degrees  of  freedom  for  both  thermal  and  electrical 
components. 

4.2.  Cylindrical  cell 

In  this  study,  we  use  the  same  battery  components  as  in  the 
first  study  rolled  into  a  cylindrical  cell.  The  cell  was  2  cm  in 
diameter  and  5  cm  long  with  eight  spiral  wounds  of  the 


Fig.  8.  Temperature  distribution  (in  I<)  in  unrolled  cell  with  consideration  of  current 
collectors:  (a)  end  of  discharge  at  1C;  (b)  2C  (Not  to  scale). 


sandwich  configuration  containing  the  same  active  materials  as 
in  the  case  of  unrolled  cell.  For  this  geometry  and  material 
composition,  the  discharge  at  1C  rate  corresponds  to  current 
density  of  22  Ah  m-2.  Material  constants  and  parameters  listed  in 
Tables  1  and  2  were  used  in  this  simulation.  Fig.  9  shows  the 
geometry  and  the  finite  element  mesh  used  to  resolve  the  ge¬ 
ometry  of  the  cylindrical  cell  including  the  current  collectors.  The 
model  has  168  zones  that  are  partitioned  into  4  quadrants.  Each 
quadrant  has  28  distinct  zones  representing  the  cell-sandwich  (S 
domain)  and  positive  (PC)  and  negative  current  (NC)  collectors. 
The  cutout  in  Fig.  9  shows  the  distribution  of  zones  within  one 
quadrant  with  electrochemical  (charge  transfer)  zones  repre¬ 
senting  cell  sandwich  shown  with  different  colors.  The  simula¬ 
tion  uses  56  concurrent  DualFoil  simulations  for  different  cell- 
sandwich  zones  within  the  cell.  Convective  (Robin)  boundary 
condition  was  applied  to  the  cell  surface  with  the  convective  heat 
transfer  coefficient  equal  to  25  W  nrr2  K_1  and  in  this  case,  the 
cell  casing  was  not  modeled.  The  results  of  temperature  distri¬ 
bution  at  the  end  of  2C  discharge  are  shown  in  Fig.  10.  The 
maximum  temperature  occurs  at  the  cell  core  as  expected.  The 
total  time  required  to  run  this  single  simulation  is  approximately 
10  min  with  a  finite  element  mesh  of  30,000  degrees  of  freedom 
for  both  thermal  and  electrical  components. 

4.3.  Pouch  cell:  experimental  validation 

In  this  study,  we  use  all  the  physics  components  for  modeling 
the  performance  characteristics  of  the  pouch  cell.  The  cell  under 
consideration  is  a  70  mm  x  110  mm  x  10  mm  4.3  Ah  pouch  cell 
manufactured  by  Farasis  Energy,  Inc.  The  experimental  discharge 
profiles  of  this  cell  for  different  C-rates  are  shown  in  Fig.  11. 

NTG  model  was  chosen  to  represent  the  electrochemical 
component  within  the  framework.  The  experimentally  obtained 
discharge  data  (Fig.  11)  was  used  to  calculate  the  variables  within 
the  linearized  polarization  model  Eq.  (1)  and  obtain  the  fitting 
constants  based  on  the  procedure  described  earlier  (Fig.  2).  Co¬ 
efficients  for  the  6-th  order  polynomial  fits  are  provided  in  Table  3. 

The  resulting  functions  for  the  open  circuit  potential  U  and  the 
polarization  parameter  Y  are  shown  in  Fig.  12.  The  coefficient  of 
determination,  R2,  for  these  fits  is  0.962  for  Y  and  0.998  for  U.  The 
simulated  cell  potential  as  a  function  of  the  depth  of  discharge  ( 6 )  is 
shown  in  Fig.  12(c)  for  1C  and  5C  applied  currents. 


Fig.  9.  Finite  element  mesh  and  distribution  of  zones  in  the  cylindrical  rolled  cell. 
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Fig.  10.  Temperature  distribution  in  a  rolled  cylindrical  cell  subject  to  2C  discharge 
current. 


The  pouch  cell  in  the  current  study  contained  17  cathode  and  17 
anode  layers,  therefore  the  finite  element  mesh  for  thermal  and 
electrical  analyses  needed  to  be  more  refined  to  resolve  the  current 
collectors  and  the  cell  sandwich.  The  finite  element  mesh  was 
divided  into  71  corresponding  zones  for  cell  sandwich  (S  domain) 
and  current  collectors  (NC  and  NP  domains).  In  addition,  two  zones 
were  added  on  the  sides  of  the  cell  to  represent  the  pouch  material 
(P  domain).  The  cell  pouch  is  coupled  to  thermal  physics  only  as  it 
does  not  participate  in  any  electrochemical  processes  or  electrical 
transport  within  the  cell.  The  thermal  conductivity  of  the  pouch 
material  was  estimated  based  on  a  three-layer  pouch  material 
structure  (DNP  Inc.,  Japan)  with  two  isolative  layers  (Nylon  and 
Polypropylene)  sandwiching  the  aluminum  layer  in  between.  Based 
on  the  parallel  and  series  formulations  for  anisotropic  thermal 
conductivity  Eq.  (10),  the  conductivity  of  the  pouch  in  transversal 
direction  was  calculated  as  kz  =  0.263  W  m-1  K_1  and  in  the  other 
two  directions  as  kx  =  ky  =  7.75  W  m_1  K-1.  Density  and  heat 


Table  3 

Coefficients  for  NTG  model  (Eq.  (1)). 


U[V] 

Y  (S  m-2) 

a0  =  4.121 

b0  =  3.6433  x  102 

ai  =  -1.0206 

bi  =  -2.0629  x  103 

a2  =  -1.1902 

b2  =  1.1248  x  104 

a3  =  8.2938 

b3  =  -2.8978  x  104 

a4  =  -20.6735 

b4  =  4.0999  x  104 

a5  =  24.8684 

b5  =  -3.1353  x  104 

a6  =  -11.0442 

b6  =  9.8589  x  103 

capacity  of  the  pouch  material  were  1150  kg  m-3  and 
1900  J  kg-1  K”1  respectively.  The  finite  element  mesh  for  the  cell  is 
shown  in  Fig.  13.  A  cutout  in  Fig.  13  represents  the  details  of  the 
distribution  of  the  zones  within  the  cell.  An  idealized  weld  with 
zero  contact  resistance  is  assumed  for  joining  all  the  current  col¬ 
lectors  to  the  tab  leads  where  a  current  flux  boundary  condition  is 
imposed.  The  simulation  uses  34  concurrent  NTG  model  instances 
for  solving  the  electrochemical  potential  in  each  cell  sandwich. 

The  results  of  temperature  distribution  at  the  end  of  two 
different  discharge  rates  are  shown  in  Fig.  15.  A  mixed  boundary 
condition  is  used  on  the  pouch  surfaces  for  the  thermal  model  to 
simulate  cooling.  Typical  values  of  heat  transfer  coefficient  corre¬ 
sponding  to  free  convection  with  air  range  from  5  to  25  W  m-2  K-1. 
The  value  for  convective  heat  transfer  coefficient  is  estimated  from 
the  temperature  measurements  of  the  cell  during  the  rest  periods 
after  each  discharge  (Fig.  14).  In  the  absence  of  the  heat  sources  and 
assuming  uniform  temperature  across  the  cell  thickness,  the  solu¬ 
tion  of  the  thermal  transport  equation  has  a  simple  exponential 
form  T  -  Tamb  =  exp (— t/r)  where  z  =  mCp/hTAs  and  hj  is  the  heat 
transfer  coefficient,  m  is  the  mass  of  the  cell,  and  As  is  the  cell 
surface  area.  The  experimental  temperature  measurements  during 
cell  rest  were  fitted  with  this  exponential  function. 

The  cooling  slope  was  very  consistent  between  the  curves  and 
the  average  was  Tavg  =  463.02  s,  which  corresponds  to 
15.48  W  m-2  I<-1  heat  transfer  coefficient.  For  the  present  nu¬ 
merical  study,  a  convective  heat  transfer  coefficient  of 
15  W  nrr2  K_1  was  used  in  the  Robin  (mixed)  boundary  conditions 
to  represent  the  natural  convection  cooling  of  the  battery  surface. 
This  configuration  results  in  a  temperature  gradient  in  the  cell’s 
through  thickness  direction  and  the  maximum  temperature  at  its 
center.  In  order  to  validate  the  modeling  results,  thermal  mea¬ 
surements  were  performed  on  the  cell  during  discharging  at 
different  C-rates  (Fig.  11).  The  temperature  was  measured  using  an 
IR  camera  and  the  data  was  obtained  at  the  geometrical  center  of 
the  cell  surface. 

A  comparison  between  simulated  and  experimentally  measured 
temperature  is  shown  in  Fig.  16.  The  markers  represent  the 
experimental  data  and  the  solid  lines  correspond  to  the  simulation 
results.  As  can  be  seen  the  predictions  of  the  temperature  values  are 
in  a  good  agreement  with  the  experiment.  The  maximum  error  of 
the  temperature  prediction  using  NTG  model  is  below  5%  for  all  of 
the  discharge  rates  considered.  The  total  time  required  to  run  this 
single  simulation  is  around  6  min  with  a  finite  element  mesh  of 
30,000  degrees  of  freedom  for  both  thermal  and  electrical 
components. 


5.  Discussion 

The  primary  goals  of  the  current  report  are:  a)  to  demonstrate 
the  flexibility  of  the  new  open  architecture  software  computational 
framework  developed  within  the  CAEBAT  program  and  b)  to  vali¬ 
date  the  simulation  results  from  the  new  3D  approach  that  resolves 
the  current  collectors  with  the  experimental  data.  The  framework  is 
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Fig.  12.  Parameters  used  in  NTG  model:  (a)  open  circuit  potential;  (b)  cell  polarization  parameter;  (c)  calculated  cell  potential  as  a  function  of  depth  of  discharge. 


designed  for  coupling  different  physics  for  thermal  and  perfor¬ 
mance  modeling  of  lithium-ion  batteries.  The  modular  structure  of 
the  framework,  where  each  physical  process  is  represented  by  a 
computer  code  component  provides  convenience  in  switching  be¬ 
tween  the  models  representing  different  processes  in  battery 
simulation.  Indeed,  in  order  to  substitute,  add,  or  remove  compo¬ 
nents  from  the  overall  model,  the  end  user  merely  needs  to  change 
a  few  lines  in  the  simulation  configuration  file. 

At  the  current  stage  three  components  have  been  integrated 
into  the  framework  representing  electrochemistry,  thermal  and 
electrical  transport.  Two-way  loose  coupling  between  the  compo¬ 
nents  was  used  in  the  present  simulations  (Fig.  5).  The  exchange  of 
the  data  between  the  components  is  set  up  in  such  a  way  (Fig.  4) 
that  the  solution  of  electrochemical  module  provides  the  re¬ 
sistances  within  the  cell  which  are  passed  to  the  electrical  module 
for  solution  for  potential  in  the  battery.  The  heat  sources  computed 
by  electrochemical  module  with  inclusion  of  the  ohmic  heat  from 
the  electrical  module  are  passed  to  the  thermal  component  to 
obtain  the  temperature  distribution.  In  this  manner,  a  true  3D 
distribution  of  potential  and  temperature  is  obtained  as  a  result  of 
the  simulation  and  as  far  as  we  know  this  is  the  first  imple¬ 
mentation  of  this  nature. 

As  in  many  instances  of  mathematical  modeling  work,  deter¬ 
mination  of  material  constants  and  system  parameters  necessary 
for  conducting  a  simulation  presents  a  challenging  task.  In  this 
regard,  the  DualFoil  model,  being  based  on  the  porous  electrode 
theory,  requires  significant  number  of  material  parameters,  many 
of  which  are  not  readily  available.  The  model  calls  for  such  pa¬ 
rameters  as  Li-ion  diffusivities  in  all  of  the  members  constituting 
the  cell  sandwich  as  well  as  composition  volume  fractions  and 
specific  surface  areas  of  the  electrodes.  The  latter  can  be  estimated 
[23]  by  assuming  the  electrode  particles  being  spherical  of  equal 
radius  Rs  as  a  =  3(1  -  e)IRs,  but  this  formula  invokes  determination 
of  electrolyte  volume  fraction  e  which  requires  porosimetry  ex¬ 
periments  being  conducted  on  the  electrode.  This  formula  was  used 
in  the  current  DualFoil  calculations  with  electrolyte  volume 


fractions  in  the  electrodes  taken  from  the  literature  (Table  1). 
Finally,  the  exponent  p  in  the  equation  for  ionic  conductivity  of 
electrolyte  Eq.  (13)  is  determined  from  the  best  fitting  the  discharge 
curve  data  [23].  NTG  model  on  the  other  hand,  while  unable  to 
provide  details  of  spatial  distribution  of  variables  (e.g.,  Li  concen¬ 
tration  in  the  cell  sandwich)  requires  only  availability  of  the 
discharge  curves  at  several  discharge  rates.  Being  empirical  in  na¬ 
ture  it  gives  good  predictions  of  the  cell  electrochemical  behavior 
and  since  all  the  parameters  are  represented  as  polynomials,  the 
savings  in  compute  time  are  tremendous.  Therefore,  this  model 
could  be  suitable  for  computations  involving  large  number  of 
concurrent  runs  of  electrochemical  component  -  that  is  for  large 
systems  such  as  battery  modules  and  packs.  The  model  requires 
determination  of  cell  conductance  (Y)  and  open  circuit  potential  ( U ) 
as  functions  of  the  depth  of  discharge.  It  should  be  re-iterated  that  Y 
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Time  /  s 

Fig.  14.  Temperature  drop  in  the  cell  during  the  rest  periods  after  discharge  at  different 
rates. 


Fig.  16.  Comparison  between  experimentally  measured  and  predicted  by  modeling 
temperatures  at  the  pouch  cell  surface. 


is  the  fundamental  electrochemical  quantity  containing  informa¬ 
tion  about  electrode  polarization  resistance  as  well  as  resistance 
of  the  electrolyte  and  separator.  This  parameter  can  be  obtained  via 
rigorous  mathematical  derivations  [35]  and  in  general  is  a  non¬ 
linear  function  of  the  cell  potential.  The  representation  in  the 
NTG  model  can  be  considered  as  a  linear  approximation  to  the 
full  polarization  curve,  derived  from  the  empirical  data.  While 
approximate,  such  treatment  yields  reasonable  predictions,  as 
can  be  seen  in  Fig.  16.  Since  the  NTG  model  describes  combined 
effect  of  processes  occurring  in  different  components  of  the  cell 
sandwich,  addition  of  the  degradation  related  models  which  are 
microstructure-dependant  (i.e.  stress  build-up  and  fracture,  SEI 
formation,  manganese  dissolution  from  LiM^CU,  etc.)  would  pre¬ 
sent  a  problem,  which  can  be  considered  as  the  major  drawback  of 
NTG  model. 

Rigorous  treatment  of  the  energy  balance  in  an  electrochemical 
cell  requires  consideration  of  heat  coming  from  polarization  (irre¬ 
versible  heat  loss),  entropy  change  (reversible  heat  loss),  as  well  as 
ohmic  losses  in  electrodes,  electrolyte  and  current  collectors  [27].  A 
detailed  treatment  of  each  of  the  components  and  discussion  of 
their  influence  on  the  combined  heat  generation  is  available  in  Ref. 
[29].  Among  the  terms  entering  the  equation  for  the  heat  genera¬ 
tion  (Eq.  (8)  and  its  simplified  form  Eq.  (9)),  determination  of 
the  reversible  heat  presents  the  major  challenge  as  it  requires 


Fig.  15.  Temperature  distribution  in  pouch  cell  at  the  end  of  discharge:  (a)  1C  applied 
current:  (b)  5C  applied  current. 


determination  of  dU/OT.  This  quantity  is  system-specific  and  is  not 
readily  available,  which  means  this  term  is  often  excluded  from  the 
modeling,  especially  when  simulations  are  performed  at  high 
discharge  rates.  The  importance  of  the  reversible  heat  obviously 
diminishes  at  higher  C-rates  when  high  overpotentials  dictate 
the  overall  heat  generation.  Elowever,  at  low  and  intermediate 
discharge  currents,  the  entropic  heat  can  constitute  as  much  as  50% 
of  the  total  heat  in  some  redox  systems  [36].  Determination  of  dUI 
dT  can  be  based  on  experimental  measurements  of  the  entropy 
change  (AS)  using  electrochemical  thermodynamic  measurement 
system  (ETMS)  technique  [37].  The  reversible  heat  source  is  related 
to  the  entropy  change  as 

-  -T§  -  74 SjL  (15) 

Where,  J  is  the  current  density,  T  is  temperature,  F  is  the  Faraday’s 
constant,  and  n  is  the  number  of  electrons  transferred  during  the 
reaction.  For  most  Li-ion  systems  n  =  1.  The  entropy  change  is 
highly  system  specific,  such  that  LiCo02  presents  much  higher 
entropy  change  than  LiFeP04  cathodes;  measurements  of  the  full 
cell  configurations  show  that  the  entropy  change  is  highly  specific 
to  the  electrode  pair  being  used  as  well  [37]. 

Naturally,  the  above  limitations  as  well  as  advantages  need  to  be 
considered  for  each  simulation  set  up  when  the  particular  models 
need  to  be  chosen  to  represent  the  physical  phenomena  in  the  best 
way.  For  instance,  as  was  mentioned,  thermal  modeling  of  large 
systems  may  be  coupled  to  NTG  model  representing  electrochem¬ 
istry  in  order  to  save  the  compute  time.  Porous  electrode  ( DualFoil ) 
model  can  be  used  when  the  details  of  distribution  of  internal 
variables  are  of  essence,  for  instance  when  coupling  with  degra¬ 
dation  models  (next  step  in  the  current  framework  development)  is 
desired.  The  structure  of  the  computational  framework  presented 
in  the  current  investigation  provides  the  user  with  convenience  of 
swapping  the  models  based  on  such  considerations. 

6.  Conclusions 

We  presented  a  new,  Open  Architecture  Software  (OAS) 
computational  framework  for  modeling  lithium-ion  batteries.  The 
framework  integrates  electrochemical,  thermal,  and  electrical 
models  for  studying  battery  performance  under  various  discharge 
rates  and  other  conditions.  Several  case  studies  representing 
different  Li-ion  cell  configurations  are  presented.  The  case  studies 


886 


S.  Allu  et  al  /  Journal  of  Power  Sources  246  (2014)  876-886 


demonstrate  the  computational  strategy  and  usage  of  different 
models  representing  electrochemical  behavior.  Importance  of  heat 
dissipation  via  metal  current  collectors  is  shown  with  the  example 
of  unrolled  cell,  where  the  temperature  can  rise  up  to  70  °C  at  2C 
discharge  when  current  collectors  are  not  included  in  the  model. 
The  validation  of  the  modeling  approach  is  presented  by  compar¬ 
ison  of  the  predicted  transient  temperatures  in  a  pouch  cell  with 
those  experimentally  measured  at  different  discharge  rates.  A  good 
correlation  between  modeling  predictions  and  experiment  was 
observed. 
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Nomenclature 


a:  specific  area 

a,:  OCP  fitting  constants  in  NTG  model 

b{.  polarization  parameter  fitting  constants  in  NTG  model 

ce:  Li-ion  concentration  in  electrolyte 

cs:  Li-ion  concentration  in  solid  phase 

Cp:  heat  capacity 

:  effective  diffusivity  in  electrolyte  phase 
Ds:  diffusivity  in  solid  phase 

/a:  mean  molar  activity  coefficient  of  the  salt  in  electrolyte 

F:  Faraday’s  constant 

h:  thickness  of  the  cell 

hj:  heat  transfer  coefficient 

i°:  exchange  current  density 

ie,  i:  current  density  in  electrolyte  and  solid  phase 

jn:  Li-ion  flux  across  the  solid— electrolyte  interface 

J:  current  density  transferred  from  the  negative  electrode  to  the  positive  electrode 

lq:  components  of  the  thermal  conductivity 

keF  effective  ionic  conductivity  of  electrolyte 

q:  volumetric  heat  source 

R:  universal  gas  constant 

Rs:  active  material  particle  radius 

r;  resistance 

t°:  lithium  transference  number 
T:  temperature 

U:  open  circuit  potential  (OCP) 

Vp:  positive  electrode  potential 
Vn:  negative  electrode  potential 
aa:  anodic  transfer  coefficient 
ac:  cathodic  transfer  coefficient 

e:  volume  fraction  of  electrolyte  in  composite  electrode 

<pe:  electrolyte  potential 

0S:  solid  phase  potential 

77;  overpotential 

p:  material  density 

op  components  of  electric  conductivity 
6:  depth  of  discharge 


